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^SJ . Abstract 

I offer a simple and useful formula for the resolvent of a small rank perturbation of 
' large matrices. I discuss applications of this formula, in particular, to analytical and 

numerical solving of difference boundary value problems. I present examples connected 
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'/^ ' with such problems for the difference Laplacian and estimate numerical efficiency of 
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the corresponding algorithms. 



^ • 1 Introduction 
in 

' Wide application of various versions of M G Krein's formula resulted lately in a marked 

. progress of the theory of boundary value problems for equations of mathematical physics. 

ps| \ In its initial form, this formula connects the resolvents of two different selfadjoint exten- 

■ sions of a given symmetric operator with finite defect indices in the Hilbert space j^]. Using 

. this connection one can find, in particular, the exact solution of the Schrodinger equation 

with the point-wise potentials, construct the correct theory of the boundary problems for 
' the Laplace operator on graphs |§, ^, obtain an approximate expression for the resolvent 

. and the exponent of the second order differential operator in the domains of M" in terms of 

the parametrix of the operator in the whole space ^] , etc. The formula for the resolvent 
of the boundary value problems, which appears in such a way, can also be used for the 
construction of numerical algorithms. 

It is interesting to find a direct analog of M G Krein's formula for the difference equa- 
' tions (and here such an analog is offered). This is important for construction of numerical 

algorithms when we reduce the initial differential equation to some finite dimensional (i.e., 
matrix) problem. If we had a discrete analog of M G Krein's formula, we would have been 
able to develop effective numerical algorithms for solution of the difference boundary value 
problems. 

I begin with a very simple matrix relation. For reasons which I will try to explain in 
what follows, I refer to it as finite dimensional analog of the Krein formula. After brief 
discussion, I present a few examples, showing the simplest applications of this relation. In 
particular, I give a short description of the algorithm for solving general boundary value 
problem for difference Laplacian in two-dimensional rectangular domain with (perhaps) 
small defects such as holes or cuts. 
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2 Low dimensional perturbations of the matrix 
and their resolvents 

We start with a very simple question from Linear Algebra. Though rather important, as 
I will try to illustrate, it did not attract attention of researchers, at least, I could not find 
it in the literature. Even Prasolov's encyclopedia of nice problems in Linear Algebra |14] 
missed it. Namely, suppose ^ is a square matrix and we had computed its inverse A~^. 
And — such a terrible but common disaster! — we observe that the typist has typed 
one of matrix elements of the initial matrix wrong! Must we redo the whole work (this is 
pretty expensive for large matrices!) or there is a cheaper possibility to obtain the correct 
answer? 

I will show an almost obvious way to answer the last half of the question in affirmative: 
there is a cheaper way. I assume, of course, that all the matrices we are going to invert 
are indeed invertible. 

Thus, suppose we have to solve the equation 

{A + B)x = f, (2.1) 

where x, f £ C", / is a known vector and x is an unknown one. If we can easily solve the 
"unperturbed" equation Ax = g for any right hand side g (this means exactly that A~^ is 
known), we can re- write ( p.l[ ) in the form 

{E + A-^B)x = A-^f. (2.2) 



We introduce new unknown vector z from the relation z = Bx. Then, if we multiply (2.2) 
by B from the left, we obtain the equation for z: 

{E + BA-^)z = BA-^f. (2.3) 

The following evident assertion holds. 

Proposition 2.1. If both matrices A and A + B are invertible, then E + BA^^ is also 
invertible. 

Proof. Assume the contrary. This means that the homogeneous equation 

{E + BA~^)z = ^ (2.4) 

has a nontrivial solution zq ^ 0. Since A is invertible by hypothesis, the vector xq = A~^zq 
exists and is nonzero. Then having substituted zq = Axq into ( |2.4| ) we obtain 

{A + B)xo = 0, 

in contradiction with the fact that A + B is invertible. I 
We can, therefore, solve equation (|2.3| ) and write 



[E + BA-^] ^BA-^f. (2.5) 
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It is clear that equation (|2.2| ) can be now rewritten as 

X = A-^f - A-^z. (2.6) 



Now, we can substitute expression ( p.5| ) for vector z into (2.6) and finally obtain the 
solution of (^.iD in the following strange form: 



x = A-^f-A-^{E + BA-^) ^BA-^f. (2.7) 
We can also rewrite this formula as an operator relation: 

~ ^ (2.8) 



[A + B)-^ = A^^ - A-^ (E + BA-^) ^ BA'^ 



which I refer in what follows as the finite dimensional analog of Krein's formula. 

Obviously, all the above is meaningless for generic matrices A and B, because the 
calculation of the inverse matrix for E + BA~^ is of the same complexity as that of the 
initial one. The situation changes dramatically if the rank of B is small. 

The word "small" means in this context that the ratio ^-r^ is much smaller than 1. In 
this case the calculation of the inverse matrix of -E + BA~^ becomes simple. 

Example 2.1. Let V be an n-dimensional vector space, e £ V and f £ V* . Let B = f 
be the linear operator in V of rank 1. Then it is easy to see that 

{E + BA-^y^ = E-{l + f {A-^e)y^ e ® {A'^ f, (2.9) 

and, therefore, we need only about operations to calculate matrix (^E + BA^^^ ^ in- 
stead of about in the general case. 

One can obtain the estimate + mn^ for complexity of such a calculation when 
rkB = m. This value is much smaller than provided ^ <C 1. So we can consider 

relation ( p.8[ ) as a version of perturbation theory in which the ratio of the rank of 
perturbation to the rank of the unperturbed operator plays the role of small parameter. 

We will see in what follows that in some important cases the calculation of the unper- 
turbed resolvent may turn out to be incredibly simple, in distinction with direct calculation 
of the perturbed one, and in these cases application of formula (|2.8| ) becomes very effective. 



Before demonstrating possible applications of formula ( |2.8D in computational mathe- 
matics, let me briefly explain the reason to baptize a very simple relation from Linear 
Algebra with a famous name. As it is mentioned in Introduction, the "actual" Krein 
formula connects the resolvents of two different self-adjoint extensions Ai and A2 of a 
given symmetric operator Aq in an infinite dimensional Hilbert space H. 

Unfortunately, it is very difficult, or, perhaps, even impossible, to read any operator 
sense into the difference Ai — A2 of such extensions, because, as a rule, this difference van- 
ishes on the intersection of their domains. Such and similar difficulties, however, had never 
been an obstacle for physicists, and they eagerly used Dirac's (5-function as a potential 
of "point- wise interaction" in the Schrodinger equation, see, e.g., [Q]. Certain arguments 
which I skip convinced me that "point-wise" perturbations of differential operators are in 
some sense perturbations of finite rank and, due to this fact, the corresponding problems 
have exact solutions. 



A Finite Dimensional Analog of the Krein Formula 



449 



Note in this connection that, although the actual Krein formula presupposes finiteness 
of defect indices of the initial symmetric operator, it may by used as well in the case of 
infinite indices (there are a number of papers on this topic, see e.g., ^] and references 
therein). Such a case arises, e.g., if we consider different boundary value problems for 
given symmetric partial differential operator 

In what follows we consider a problem of calculating resolvents of extensions (see formal 
definition in the next section) of difference operators based on formula (^.§|). 

It seems that the corresponding relation is a finite dimensional analog of the relation 
for differential operators. 

This impression is not an illusion. Indeed, it is possible to consider (in some well- 
defined sense) the difference operators as approximation of differential ones, and then 
one can prove that in the case of finite defects (e.g., for ordinary differential operators) 
our formulas converge to the corresponding formulas for differential operators (private 
communication of E Gordon and S Albeverio; together with them we intend to explain 
this in detail elsewhere). Such convergence plays a crucial role both for goals of numer- 
ical analysis and as an instrument for investigation of infinite dimensional operators via 
their finite dimensional approximations (see, e.g., |5|). Unfortunately, rigorous results 
about convergence of finite-dimensional approximations of operators requires for proofs 
a nonelementary technique which is out of frame of this work. Nevertheless, I consider 
(briefly and without proof) at the end of Section ^ the simplest example of convergence of 
finite dimensional Krein formula for difference approximations of the Schrodinger operator 
with 5-potential on the unit circle to the usual Krein's formula for the resolvent of this 
operator. 

3 Boundary-value problem for the difference operators 

Difference approximations of boundary value problems for differential operators are a 
base for the numerical solving of such problems. Here I just introduce a convenient for 
our nearest goals language for formal description of "abstract" difference boundary value 
problems. I could not find an appropriate analog of such a language in the literature. 
Hopefully, the following examples make it clear why this language is useful and convenient. 

For any set M, let C{A4) be the space of all complex- valued functions on M. 

A linear map A : C{M) — > C{M) will be called a formal difference operator in C{M) 
if for each x £ M. there exist a finite set "yA{x) C Ai and function £ (7(7^ (x)) such 
that 

{Af){x)= f^c{M), (3.1) 

ye7A(a;) 

Let be a subset of M.. The point 2; G is an inner point of the set Q with respect to 
the map A, if 7a(-z) ^ ^• 

The point 2: E $7 is a boundary point of the set Q with respect to the map A, if 7^1(2:) \r2 ^ 
0. The boundary of 0, with respect to A is the set of all boundary points of 0,. Define 
the set of exterior points of 0, with respect to A to be 

bAn= \J (7^(x)nn), where Ti = M\n. 
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Note that in the "difference" case the sets and 6^0 do indeed depend on the map A 
in contrast with the continuous situation. Observe that 

y 7a(x) = nu bAn. 

Define the map Aq : C(i7 U Ia^) — > C(i7) by formula ( |3.lD for any point x E ^2. 

Let L : C{0,) — > C(riU6^il) be a hnear map such that {Lf)[x) = f{x) for ah x G O. 
The operator L will be called an extension operator for the map A. 

We say that the operator A^ : C{Q) — > C(il) is an L-extension of the formal difference 
operator A if 

{ALf){x) = {AnLf){x), xen. (3.2) 

Note, that the in case described, the extension operators play the role of boundary 
conditions for differential operators. I hope that this will be clear from the examples of 
this section. 

In what follows I suppose that the set 0, is finite. Let for A € C the operator -Rl(A) be 
the resolvent of the L-extension of the operator A such that RlW = (A^ — XE)~^ . 

We show that formula (^]^) establishes a simple algebraic connection between resol- 
vents RlW and Rk{X) of two different extensions of the formal difference operator A 
corresponding to two extension operators L and K. (In what follows I assume that A is 
a common resolvent point for both Al and Ak-) To obtain such a connection, note first 
that definition ( |3.2| ) implies 

Ak = Al + Dlk, where Dlk = An{K - L). (3.3) 

Now, let us replace matrix A^^ in (|2.8|) with Rl{X), matrix (^4 + B)~^ with Rk{X), and 
B with Dlk- Then we see that 

Rk{X) = Rl{X) - Rl{X){E + DLKRL{X)r^DLKRL{X). (3.4) 

Observe that all the inverse operators in this formula exist by the hypothesis. 
What do we gain from this formula? Note first of all, that it is easy to see that 

rkDLK <#{dA^), and vk{AK - X) = vk{AL - XE) = #{n). 

It is remarkable that as a rule (see examples in what follows) ^{Oa^) *C #(f^), and we 
are in the situation discussed in Section |l|. Examples also show that the complexity of 
calculation of the resolvent for different extensions may be essentially different. 

Example 3.1. Resolvent of the one-dimensional difference Laplacian. Although 
it seems that this example has no practical meaning, it makes very clear all previous 
abstract constructions and has all essential features of practically important Example 3.2. 

Let = Z the set of integers, a positive integer, and O = {0, 1, . . . , A^— 1}. The one- 
dimensional difference Laplacian is the formal difference operator A : C(A4) — > C{M) 
defined by the relation 



(A/)(x) =/(x + l)-2/(x) + /(x-l), where xeM and f e C{M). 
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Note that operator A differs from the usual difference approximation of the differential 
expression on the uniform grid in M} by a factor only. 

In the case considered d/^Vl = {0, — 1} and b/^Vt = { — 1, N}. 

We list all the extensions for A. Let / : C(f2) — > C({— 1, A}), i.e., / is represented by 
a 2 X A complex matrix. Then for / G C(Q) set 

{Lf){x) = f{x), xen, {Lf){y) = {lf){y), y e b^n. 

Clearly, any extension operator for A must be of such form. 

Among all L-extensions of A there exists an exceptional one, for which the correspond- 
ing resolvent has an "almost explicit" expression. This is the so-called periodic extension, 
defined by extension operator Lq such that 

(iof) (-1) = /(A - 1), (fo/) (A) = /(O). (3.5) 

The exceptional role of this extension (denoted in what follows by Aq instead of A/,,, for 
brevity) is the consequence of the fact that it can be diagonalized by Discrete Fourier 
Transformation (DFT), i.e., 

Ao = FAF*, (3.6) 
where A is the multiplication operator (i.e., the diagonal matrix) 

(A/)(x) = -4sin2^/(x), xGn, 

and the unitary DFT operator F is defined by the relation 



^ yen 

Formula ( |3.6| ) immediately implies the equality 

Rq{X) = F{A- XE)-^F*, (3.7) 

and this is what we meant under the explicit formula for resolvent. 

It is well known that there exists an abnormally effective numerical method (called Fast 
Fourier Transformation, or FFT) for application of DFT to the vector. It requires only 
r-^ NlogN arithmetic operations instead of ^ for the general N x N matrices ||lOt| . 
This fact crucially reduces the complexity of computation of operator (|3.7|). 

Is there an algorithm which allows one to calculate the resolvent of an arbitrary exten- 
sion of A with the same complexity as for Aq? Formula (^^) gives a positive answer to 
this question. It only suffices to show that the computation of matrix {E + Dqk Ro{X))^^ is 
not a problem. Indeed, due to the fact that {DoKf){x) = for / G (7(0) and x £ Q\dA^, 
to solve the equation 

{E + DoKRo{X))f = 9, (3.8) 

we only have to find /(O) and /(A — 1). We denote by 6x the function from C(0) defined 
by 

II, tor X = y. 
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Since f{x) = g{x) for x ^Vt \ d/\Q, we can re-write equation (3^) in the form 



(1 + (Do/<^o(A)5o)(0))/(0) + (DoKRoWSN-imfiN - 1) 
= 5(0) - {DoKRoWig - 9{0)6o - 9{N - l)5jv-i))(0), 
{DokRo{X)So){N - l)/(0) + (1 + {DokRo(.X)Sn-i){N - l))f{N - 1) 
= g{N - 1) - {DokRo{X){9 - 9{0)So - 9{N - 1)6n-i)){N - 1). 

This is a system of two hnear equations for two unknowns, which is solvable due to 
Proposition 2.1. So, to calculate the resolvent Rk{X), we only have to know how to 
calculate Ro{X) and how to invert 2 x 2-matrices . . . 

This example is, as have already been said, of no practical importance, because there 
exists another (not DFT-based) algorithm for inverting the general three-diagonal matrix 



of complexity ~ N (so called sweep method, see, e.g., |1C]). For most often used types of 
boundary conditions (i.e., extension operators K), the matrix of Ak is of this kind, and 
the sweep method becomes preferable. For example, the Dirichlet problem corresponds to 
the extension defined by the map I of the form 

{if) (-1) = -/(o), (//) (iV) = -f{N - 1). 

and this leads to a three-diagonal matrix. 

Note, however, that if for an extension operator K the matrix of Ak is not three- 
diagonal (as is the case, e.g., for Aq), one can use the Dirichlet extension as the "initial" 
one and solve the problem for the ii'-extension using only ~ N arithmetic operations! The 
reason for using DFT in this example becomes clear from the following example. 

Example 3.2. The boundary value problem of third kind for Laplacian in two- 
dimensional rectangle. We consider now the boundary value problems for the two- 
dimensional difference Laplacian. Let Ad = Z^, N and M positive integers, and Q = 
{0, . . . , — 1} X {0, . . . , M — 1}. The formal two-dimensional Laplace operator which we 
denote by the same symbol A : C{Ai) — > C{M) is given by the formula: 

(A/)(x, y) = fix + 1, y) + f{x - 1, y) + f{x, y + l) + f{x, y - 1) - 4/(x, y), 
for x,yeM and feC{M). 

Clearly, the set of the boundary points with respect to the operator A is 

d^n = ({0, iV - 1} X {0, . . . , M - 1}) U ({0, . . . , - 1} X {0, M - 1}) , 

so that #{dA^) = 2{N + M - 2). We see once more that #{dA^) < #{^) = MN. 
Hence, there exists a good chance for applying Krein's formula. To actually apply it, we 
first describe the set 6a ^• 

The next geometric proposition is almost evident and we omit proof. 

Proposition 3.1. 1) For p = {x,y) G 1? set \p\ = \x\ + |y|. Then p G b^.^ if and only if 
there exists (and then it is unique) e{p) G 1? such that |e(p)| = 1 and p + e{p) £ Cl. 
2) #{bA^) = 2A + 2M. 
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We will not describe all extension operators for A (though possible, this is not inter- 
esting), instead we will consider several distinguished cases. First of all, as in the one- 
dimensional case considered in Example 3.1, there exists a remarkable periodic extension 
defined by extension operator L of the form 

{Lf){x,y) = f{x mod N, y mod M), (a;,y) G 6a^- 



The corresponding operator will be denoted again by Aq and it has the same characteristic 



property, namely, may be diagonalized by a two- dimensional DFT [10|. Therefore, one 
needs ^ MNlogMN arithmetic operations for calculating the resolvent Rq{X) instead of 
about (MN)^ to invert the general linear operator in C{Q). 

Among other extensions of two-dimensional difference Laplacian, I consider only the 
ones corresponding to local boundary conditions for the differential Laplace operator. These 
extensions are defined by the family of extension operators K of the form 



{Kf){p) = k{p)fip + e{p)), p E bA^, (3.9) 

where k S C{b/\i}) and e{p) is defined in Proposition 3.1. 

We now consider again the relation (^]^). It is easy to see that, as in the one-dimensional 
case, this equation can be transformed to a linear equation for function / G C((9a^) and 
we need ~ (2M -|- 2N — 2)^ arithmetic operations to solve it. For M, large enough, 
the inequality (2M + 2N - 2f < {MNf holds, and we obtain the algorithm for solving 
the third kind boundary value problem for two-dimensional difference Laplacian with 
complexity ~ (2M -|- 2A^ — 2)^ -|- MN log MN arithmetic operations. Moreover, if one has 
to repeatedly solve this problem for different right hand sides, it suffices to calculate matrix 
{E + DokRoW)"-^ only once and then we need only ~ (2M + 2N - 2f + MNlogMN 
arithmetic operations for each right hand side. Asymptotically, this complexity is the 
same as that for the periodic Laplacian. 

Note that in contrast with the one-dimensional case, the direct (i.e., non-iterational) 
methods for calculation of Rk{X) exist only for exceptional extension operators even from 



family (^.S|), see |10|. It makes Example 3.2 important in practical applications. 



Example 3.3. The Laplacian in the two-dimensional rectangle w^ith a hole. Let 

M, A and be the same as in Example 3.2 and p an inner point of Q. Let Qp = Q \ {p}. 
It is clear that p £ b^Qp and 6a = {p} U 6a^- We consider the extension operator Kp 
of the form (|3.9| ) and suppose in addition that 

{Kpf){p)= c,{e)fip + e), a{e)eC. 

{e:\e\=l} 

Note that in this formula p + e £ Qp ior all e due to our hypotheses. 

It is easy to see that operator is exactly a rank 1 perturbation of Ax|c(np)) where 
we consider the space C{^lp) as a subspace in C(0) consisting of functions / such that 
f{p) = 0. So the resolvent of A^p can be calculated with the same efficiency as that of 
Aj<-! This is indeed remarkable, because one can consider operator A^^ as the difference 
approximation of the differential Schrodinger operator with point-wise potential [Q], and 
we see that the difference case can be investigated with the help of the introduced finite 
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dimensional analog of the Krein formula in the same manner as differential operators with 
point-wise potentials are investigated by means of the "actual" Krein formula. 

It is clear that in the same way one can construct resolvents for Laplacian in rectangle 
with more complicated defects (like holes containing more than one point, cuts, etc). Our 
approach is efficient provided i^(J)A^) <C and we know an effective algorithm for 

calculating resolvent of at least one extension. 

Example 3.4. The point-wise potentials in one-dimensional case and conver- 
gence. The aim of this example is to demonstrate that in simplest case application of 
formula ( |2.8| ) to the difference approximation of differential operator leads to the expres- 
sion for resolvent which term-by-term converges to one obtained by applying the "actual" 
Krein's formula to initial differential operator. 

Let A be the Laplace operator — ^ in L2([0,27r]) with periodic boundary conditions. 
It is evident that its resolvent is of form 

1 °° 1 

{R^{X)^){x) = - J2 V'me^™^^-, (3.10) 

m=—oo 



where 



f-2n 

ipm= ip{x)e-'"'''dx. (3.11) 
Jo 



Following Krein, consider the one-parametric family of self-adjoint extensions of the 
restriction of A onto the space of smooth functions vanishes in the neighborhood of the 
endpoints of the interval [0, 2tt], such that the resolvents of operators from the family are 
of the form 




(3.12) 



where is a parameter of family. (Note that all series in this expression converge either 
in L2([0, 27r]) or in C when ip € L2([0, 27r]).) It is well-known (see, e.g., where a 
number of similar examples are considered) , that for each real fi the operator Ra^ is indeed 
the resolvent of a self-adjoint operator A^ in L2([0, 27r]). This A^ is usually called the 
Schrodinger operator with 5 -potential (parameter /x plays the role of a coupling constant)]^. 

Let now M be a positive integer and let operator Am in L2({0, . . . , 2M — 1}]) be of the 
form 

(AM/)i = -;^(/i~i-2/, + /,+i), j G {0,...,2M-1}, (3.13) 



^Of course, formula (3.12) does not give all possible extension of symmetric operator considered, but 



the family described suffices for our goals. 
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where h = vr/M and the "exterior" values of / are defined by "periodic boundary condi- 
tions" ( |3.5D . It is easy to see that the resolvent of operator Am is, due to relation (3.7), 
of the form 

M 



27r ^ 

m=-M+l 



ihmj 



hm 



h? 



sm 



where 



2M-1 

h ipji 

m=0 



-ihmj 



(3.14) 



(3.15) 



Observe that relations (|1^), (|3ll|) and (|3l|), (ISl^ ) are of similar form. Moreover, 
setting Th : C([0, 27r]) ^ L2({0, . . . , 2M - 1}) as 



one can see that that the relations 



iG{0,...,2M-l}, 



(3.16) 



hm \\ThAf-AMThf\\h = Q, 



lim \\ThRA{X)f -RAMWThf\\h = (3.17) 



Af-»oo 



hold for every sufficiently smooth periodic function / G L2{[0, 27r]) if the norm in L2({0, . . . , 
2M - 1}) is 

2M-1 

This means exactly that the family of finite dimensional operators Am, M G N, approxi- 
mates the operator A []To| , |5|, or, in another words. Am tends to A when M — > oo. 
Let now be the operator in L2{{0, . . . , 2M — 1}) given by the formula 

1 



h' 



(3.18) 



It is easy to make use of (^j^) in order to calculate the resolvent of Am^ = Am + /^^g^ 
(cf. also with Example 2.1 and relation ( p.9[ )). In this way we obtain an expression for 
the resolvent of A Mfj. ■ 

1 ^ 

m=-J\/+l 



1 

, hm 



' ^ M 

2^ ^ 



4 - 2 hm \ 

M+1 W^^"^ 




sm 



A 



^ihmj 



4 - 2 hm \ 
M+1 Z?*^™ 



1 + ^ 

27r 



A/ 

E 



(3,19) 



_ 4 - 2 hm \ 
M+1 /I^ ^™ — - 



We compare now relations (|3l^ ) and ( pO] ). It is easy to see that for resolvents and 
Ram^ a relation like ( ^.17 ) holds. Moreover, one can see also that each term in the left 
hand side of (3.19) converges to the corresponding term in (3.12). Hence, one can assert 
that in this sense the finite dimensional Krein formula converges to the "natural" Krein 
formula for the resolvent Ra ■ This fact gives an additional argument in favor of the name 
"finite dimensional analog of Krein formula" for relation (|27 
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4 Concluding remarks 

The finite dimensional analog of Krein's formula proved to be a useful instrument for 
investigation of difference equations both analytically and numerically. Moreover, it gives 
us a new approach to study differential problems (and, more generally, other "continuous" 
extensions) by reducing them to the corresponding difference (or, more generally, other 
finite dimensional) approximations. In this connection it is interesting that, in contrast 
with the "actual" Krein formula, our algebraic relation does not require operators involved 
to be Hermitian. 

The method proposed for solving difference boundary value problem is applicable to 
a wide class of equations, in particular, in the case of complicated multi-point boundary 
conditions for one-dimensional equations, for the rectangular two-dimensional domains 
with cuts and some other "small" defects, for some cases of variable coefficient of difference 
operators, etc. It is clear that in every specific case one needs to adapt the general 
algorithm described in Sections 2, 3, but this general scheme is, nevertheless, useful for 
construction of particular numerical procedures. 

The described method for solving of boundary value problems was successfully used 



in [|11|. 



A similar approach exists also for constructing other than resolvent functions of differ- 



ence operators. This is needed, e.g., in initial boundary value problems, see [12|. Certain 
moments of the method proposed for solving difference boundary value problems were 
announced in p^ . 
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